For this exercise, I chose to replicate an interactive map: https://projects.fivethirtyeight.com/redistricting-maps/

drawing

Data used in The Atlas Of Redistricting: https://github.com/fivethirtyeight/redistricting-atlas-data

Weirdly, though knit worked without errors, building the website didn't work until I substituted select() with subset() when using pipes with spatial data

#Load packages
library(raster)
library(rgdal) 
library(sf)#st_read to read shapefiles
library(maptools)
library(broom)
library(ggplot2)
library(tidyverse)
library(leaflet) #interactive map
library(tigris)
library(grid)
library(tmap) 
library(RColorBrewer)

Why are some states missing "current" district shapefiles?

According to the readme in GitHub "In cases where the current congressional map fulfilled the goals of a custom map for a state, a shapefile is not included." What shapefiles did they use then?

MAPTYPE = Current : "The current congressional boundaries"

#Load the shapefiles
baseDir = "C:/Users/Gabriella Veytsel/Desktop/MADA 2021/Gabriella_Veytsel-MADA-portfolio/data/redistricting-atlas-data-master/shp"

#There are multiple shapefiles for each state in the shp folder
#For this map, looks like they the shapefile with MAPTYPE = current 
filenames <- dir(baseDir, "*current.shp") 
filepaths <- paste(baseDir, filenames, sep='/')

listOfShp <- lapply(filepaths, st_read) #read each file into a list
all_shapefiles <- do.call(rbind, listOfShp) #combine shapefiles into one

#Load the metadata
districts <- read_csv("C:/Users/Gabriella Veytsel/Desktop/MADA 2021/Gabriella_Veytsel-MADA-portfolio/data/redistricting-atlas-data-master/districts.csv")

Found missing state shapefiles from census website, will use them to supplement the missing ones https://www2.census.gov/geo/tiger/PREVGENZ/cd/cd103shp/

baseDir_state <- "C:/Users/Gabriella Veytsel/Desktop/MADA 2021/Gabriella_Veytsel-MADA-portfolio/data/redistricting-atlas-data-master/shp_gv_census/Missing Shapefiles"

filenames_state <- dir(baseDir_state, "*.shp") 
filepaths_state <- paste(baseDir_state, filenames_state, sep='/')

listOfShp_state <- lapply(filepaths_state, st_read) #read each file into a list
all_shapefiles_state <- do.call(bind_rows, listOfShp_state) #combine shapefiles into one

#These states only have 1 congressional district
all_shapefiles_state <- subset(all_shapefiles_state, select =ST)

#Missing state names
table(all_shapefiles_state$ST)
all_shapefiles_state <- all_shapefiles_state %>%
  mutate(STATE = case_when(
         (ST == "02") ~ "AK",
         (ST == "56") ~ "WY",
         (ST == "30") ~ "MT",
         (ST == "38") ~ "ND",
         (ST == "46") ~ "SD",
         (ST == "50") ~ "VT"))

Merge metadata with each shapefile (different criteria for keys)

all_shapefiles_state_merge <- merge(all_shapefiles_state, districts, by.x = "STATE", by.y = "state") %>%
  filter(maptype == "current")

districts_merge <- merge(all_shapefiles, districts, by.x=c("MAPTYPE", "STATE", "DISTRICT"), by.y=c("maptype", "state", "district")) 

Merge shapefiles together

districts_merge <- bind_rows(all_shapefiles_state_merge, districts_merge) 
districts_merge <- subset(districts_merge, select = -c(maptype, ST, district, statefp))

Explore the variable dem_chance for choropleth

dem_chance: The probabilities of electing a Democrat or Republican are based on how often seats with a given Cook PVI elected members of each party between 2006 and 2016.

glimpse(districts_merge)
## Rows: 511
## Columns: 20
## $ STATE                <chr> "AK", "AK", "AK", "AK", "AK", "AK", "AK", "AK", "~
## $ population           <dbl> 710231, 710231, 710231, 710231, 710231, 710231, 7~
## $ population_18_over   <dbl> 522853, 522853, 522853, 522853, 522853, 522853, 5~
## $ PVI                  <dbl> -9.39, -9.39, -9.39, -9.39, -9.39, -9.39, -9.39, ~
## $ dem_chance           <dbl> 5.409427, 5.409427, 5.409427, 5.409427, 5.409427,~
## $ `Non-Hispanic White` <dbl> 68.27674, 68.27674, 68.27674, 68.27674, 68.27674,~
## $ `African-American`   <dbl> 3.084806, 3.084806, 3.084806, 3.084806, 3.084806,~
## $ `Hispanic/Latino`    <dbl> 4.67378, 4.67378, 4.67378, 4.67378, 4.67378, 4.67~
## $ Asian                <dbl> 5.332856, 5.332856, 5.332856, 5.332856, 5.332856,~
## $ `Native American`    <dbl> 13.27008, 13.27008, 13.27008, 13.27008, 13.27008,~
## $ `Pacific Islander`   <dbl> 0.8598975, 0.8598975, 0.8598975, 0.8598975, 0.859~
## $ Other                <chr> "4.501838948997137%", "4.501838948997137%", "4.50~
## $ race_category        <chr> "Non-Hispanic White Majority", "Non-Hispanic Whit~
## $ minority_chance      <dbl> 9.889405, 9.889405, 9.889405, 9.889405, 9.889405,~
## $ current_map          <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N~
## $ impossible           <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N~
## $ MAPTYPE              <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N~
## $ DISTRICT             <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N~
## $ STATEFP              <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N~
## $ geometry             <POLYGON> POLYGON ((-143.7157 70.0248..., POLYGON ((-16~
summary(districts_merge$dem_chance)
##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
##   0.00374   4.10570  10.79572  39.78434  93.75583  99.99987

In a separate article, I found information on the bottom legend of their map:

PVI (Cook Partisan Voting Index) measures how much more Democratic or Republican a district voted relative to the national result in an average of the last two presidential elections. We categorized partisanship into three "buckets":

https://fivethirtyeight.com/features/we-drew-2568-congressional-districts-by-hand-heres-how/

districts_PVI <- districts %>%
  filter(maptype == "current") %>%
  mutate(district_pvi = 
         ifelse(between(PVI, -5, 5), "HIGHLY COMPETITIVE DISTRICTS",
         ifelse(PVI < -5, "USUALLY REPUBLICAN DISTRICTS",
         ifelse(PVI > 5, "USUALLY DEMOCRATIC DISTRICTS", 0))))
table(districts_PVI$district_pvi, useNA = "always")
## 
## HIGHLY COMPETITIVE DISTRICTS USUALLY DEMOCRATIC DISTRICTS 
##                           72                          168 
## USUALLY REPUBLICAN DISTRICTS                         <NA> 
##                          195                            0
state_pvi <- districts_PVI %>%
  subset(select = c(state, statefp, district_pvi)) %>%
  group_by(state, statefp, district_pvi) %>%
  summarize(n = n())

state_pvi <- spread(state_pvi, district_pvi, n) #have to tidy dataset for all 3 rows to appear on map
state_pvi[is.na(state_pvi)] <- 0 #if missing, make 0

state_map <- st_read("C:/Users/Gabriella Veytsel/Desktop/MADA 2021/Gabriella_Veytsel-MADA-portfolio/data/redistricting-atlas-data-master/shp_gv_census/gz_2010_us_040_00_500k.shp")

states_merge <- merge(state_map, state_pvi, by.x = "STATE", by.y = "statefp")

Plot choropleth using leaflet package
Struggling to move Hawaii, tried using the grid package and st_transform
Hover over each state with mouse

pal <- colorNumeric(c("red", "royalblue"), domain=districts_merge$dem_chance) #setting up the colors
labels <- paste(
    states_merge$NAME,"<br/>", 
    "USUALLY DEMOCRATIC DISTRICTS: ", states_merge$`USUALLY DEMOCRATIC DISTRICTS`, "<br/>", 
    "HIGHLY COMPETITIVE DISTRICTS: ", states_merge$`HIGHLY COMPETITIVE DISTRICTS`, "<br/>",
    "USUALLY REPUBLICAN DISTRICTS: ", states_merge$`USUALLY REPUBLICAN DISTRICTS`,
    sep="")  %>% 
  lapply(htmltools::HTML)
  
m <- leaflet() %>%
  addProviderTiles("CartoDB.Positron") %>%
  setView(-98.483330, 38.712046, zoom = 4)

#Color by district
m1 <- m %>%
  addPolygons(data = districts_merge, 
              fillColor = ~pal(districts_merge$dem_chance), 
              fillOpacity = 0.7, 
              color = "white",
              weight = 0.2, 
              smoothFactor = 0.2) 

#Add highlight when hovering over the state, not district
m2 <- m1 %>% addPolygons(data = states_merge,
                         fillColor = "transparent",
                         weight = 0.5,
                         opacity = 1,
                         color = "white",
                         dashArray = "0",
                         fillOpacity = 0.7,
                         highlight = highlightOptions(
                           weight = 5,
                           color = "#666",
                           dashArray = "",
                           fillOpacity = 0.7,
                           bringToFront = TRUE),
                         label = labels,
                         labelOptions = labelOptions(
                           style = list("font-weight" = "normal", padding = "3px 8px"),
                           textsize = "15px",
                           direction = "auto"))
#Add legend
m3 <- m2 %>% addLegend(pal = pal, 
                       values = districts_merge$dem_chance, 
                       position = "bottomright", 
                       title = "Chance of being represented by Democratic party") 

m3